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We present a new nonempirical density functional generalized gradient approximation (GGA) that 
gives significant improvements for lattice constants, crystal structures, and metal surface energies 
over the most popular Perdew-Burke-Ernzerhof (PBE) GGA. The new functional is based on a 
diffuse radial cutoff for the exchange-hole in real space, and the analytic gradient expansion of 
the exchange energy for small gradients. There are no adjustable parameters, the constraining 
conditions of PBE are maintained, and the functional is easily implemented in existing codes. 
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Kohn-Sham density functional theory (DFT) 0, 0] 
makes it possible to solve many-electron ground-state 
problems efficiently and accurately. The DFT is exact if 
the exchange-correlation (XC) energy Exc were known 
exactly, but there is no tractable exact expressions of 
Exc m terms of electron density. Numerous attempts 
have been made to approximate Exc, starting with the 
local (spin) density (LSD) approximation (LDA), which 
is still widely used. The generalized gradient approxi- 
mations (GGAs) d- E3- 01 are semilocal, seeking to im- 
prove upon LSD. Other more complicated approxima- 
tions are often orbital-dependent or/and nonlocal. They 
suffer from computational inefficiency; it is much harder 
to treat them self-consistently and to calculate energy 
derivative quantities. 

The XC energy of LSD and GGAs are 
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is the XC energy density for the uniform electron 
gas. LSD is the simplest approximation, constructed 
from uniform electron gas, and very successful for solids, 
where the valence electron densities vary relatively more 
slowly than in molecules and atoms, for which GGAs 
[HE| achieved a great improvement over LSD. It is well 
known that LSD underestimates the equilibrium lattice 
constant a by 1-3%, and some properties such as fer- 
roelectricity are extremely sensitive to volume. When 
calculated at the LSD volume, the ferroelectric instabil- 
ity is severely underestimated 0, S 0- On the other 
hand, GGAs tend to expand lattice constants. They well 
predict correct ao for simple metals, such as Na and K 
|6j, however for other materials t hey often overcorrect 
LSD by predicting ao 1-2% bigger than experiment. 
Predicting lattice constants more accurately than LSD 
remains a tough issue, even for state-of-the-art meta- 
GGAs; nonempirical TPSS only achieves moderate 



TABLE I: Calculated equilibrium volume Vb (A 3 ) and strain 
(%) of tetragonal PbTiOs for various GGAs comparing with 
experimental data at low temperature pM|. 





PW91 


PBE 


revPBE 


RPBE 


Expt. 


Vb 


70.78 


70.54 


74.01 


75.47 


63.09 


strain 


24.2 


23.9 


28.6 


30.1 


7.1 



improvement over PBE, while empirical PKZB 12] is 
worse than PBE. GGAs are especially poor for ferro- 
electrics, e.g., PBE |5( predicts the volume and strain 
of relaxed tetragonal PbTiOs more than 10% and 200% 
too large, respectively [lj|, and other GGAs [1 Hi. 
are even worse, as seen in Table [i] Another more compli- 
cated functional, the nonlocal weighted density approx- 
imation (WDA) is also unsatisfactory for this case [l3| . 
To compute these properties correctly, people often con- 
strain volumes at their experimental values V^ xp t. How- 
ever, V Gxp t is not available for predicting new materials, 
certain properties are still wrong even at Ve Xp t; theoret- 
ically it is more satisfactory to do calculations without 
any experimental data adjustments. 

In order to study finite temperature properties, e.g., 
ferroelectric phase transitions, effective Hamiltonian and 
potential models, which are used in molecular dynamics 
(MD) or Monte Carlo (MC) simulations, have been de- 
veloped with parameters fitted to first-principles results. 
When the model parameters are fitted to LSD data, these 
simulations greatly underestimate the phase transition 
temperatures T c at ambient pressure [lfil E3- and 
overestimate if fitted to the GGA results |l7j. A sim- 
ple but more accurate approximation for XC energy is 
necessary. 

Because the magnitude of the exchange energy is much 
bigger than correlation in most cases, and PBE generally 
describes the correlation energy with enough accuracy, 
we focus only on the exchange in this paper. The dimen- 
sionless reduced gradient s = | Vn|/[2(37r 2 ) 1 / 3 n 4 / 3 ]. The 
exchange enhancement factor Fx is defined as 
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FIG. 1: (Color online) Exchange enhancement factors Fx 
as functions of the reduced gradient s. The green curve is 
Fx of TPSS meta-GGA with a = 1.0, which corresponds 
to slowly varying density limit. Symbols in the inset are Fx 
determined by the real-space cutoff procedure. A sharp cutoff 
generates the crosses, and its parameters are chosen to fit 
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The circles correspond to a diffuse cutoff with the 



same parameters. 



where the exchange energy density of the uniform elec- 
tron gas e x nil (n) = -^-(Z^n) 1 / 3 . The PBE ensatz of 
Fx has the general form 



F x = 1 + k- k/(1 + x/k), 



(4) 



where k = 0.804 to ensure the Lieb-Oxford bound [l9j, 
and x = /is 2 with p, = 0.21951 to recover the LSD linear 
response, i.e., as s — > 0, the exchange gradient correction 
cancels that for correlation. In the range of interest for 
real systems < s < 3, the PBE Fx is a simple numerical 
fit to that of PW91, which is constructed from the gradi- 
ent expansion of a sharp real space cutoff of the exchange 
hole 0,113, plus some exact constraints ||. Unlike atoms 
and molecules, solids can have a diffuse tail around the 
exchange-correlation hole, and a diffuse radial cutoff fac- 
tor [1 + (u/u x ) 2 ] exp[— (u/u x ) 2 ), where u is the distance 
from the hole center and u x is the fixed radial cutoff, leads 
to a smaller Fx for s > 1 than that of the sharp ra- 
dial cutoff (inset of Figure^) . This explains why the PBE 
(PW91) functional improves total energies of atoms and 
atomization energies of molecules greatly over LSD, but 
often overcorrects LSD for solids. Two revised versions 
of PBE, namely revPBE EH with empirical k = 1.245 
and RPBE with F x = 1 + k - k exp(- fis 2 / k) , fur- 
ther exaggerate Fx (Figure P), giving better energies for 
atoms and molecules, but worse lattice constants of solids 
(Table HJ). 

The real space cutoff procedure can only give quali- 
tative features of Fx, not the exact behavior because it 
depends on the detailed approximations of the procedure 
and fitting to parameters, so other known constraints 



must be chosen to determine Fx- Usually valence elec- 
tron densities of solids vary much more slowly than elec- 
tron densities of atoms and molecules. The choice of p, in 
PBE violates the known gradient expansion of Svendsen 
and von Barth |2lT | for slowly varying density systems, 

where p = s 2 , q = V 2 n/[4(37r 2 ) 2 / 3 n 5 / 3 ] is the second or- 
der reduced gradient, and D = is the best numerical 
estimate. If p is set to |j, Fx of equation 0] will be low- 
ered. However the behavior of PBE Fx for small s needs 
to be retained because (i) it is necessary to retain cancel- 
lation of gradient correction of exchange and correlation 
as s — y 0; (ii) Fx determined by a diffuse radial cutoff is 
close to that by sharp radial cutoff for small s (inset of 
Figure Thus we propose the following ensatz for x in 
Equation (@J: 
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where the parameter c is set to recover the fourth order 
parameters in Equation (JSJ) for small s. Because a good 
ation of q for slowly varying densities is q « |p 
146 (|) 2 -^I + (M-^) = 0.0079325. Our new 
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functional will be referred to as "WC" , and F^ still sat- 
isfies the four conditions (d)-(g) constraining F^ BE [|J- 
It is nearly identical to Fx BE for s < 0.5, and smaller 
for bigger s, as displayed in Figure ^ Amazingly, the 
present simple Fx C matches that of the more byzantine 
TPSS meta-GGA for slowly varying densities very well, 
which makes use of the kinetic energy density to enforce 
Equation JSJ for small p and q. Contrasting with earlier 
attempts such as RevPBE, WC has no adjustable param- 
eters fitted to experimental data, and it is constructed 
completely from fundamental physical considerations. 

We tested the new functional of Equation © by com- 
puting equilibrium crystal structures and cohesive ener- 
gies of solids, jellium surface energies, and exchange en- 
ergies of atoms. We used the planewave pseudopotential 
method (ABINIT4.4.4 H2) for solids, and an all-electron 
atomic code for atoms. It is straightforward to imple- 
ment the current GGA from the PBE pseudopotential 
code. We used the same configurations for each atom to 
generate optimized norm-conserving pseudopotentials by 
the OPIUM code of LSD, PBE, and WC. 

First we calculated equilibrium lattice constants ao of 
18 solids as tested in Refs. [H Ej- We found LSD 
underestimates, while PBE overestimates ao, and current 
LSD and PBE errors agree well with previous ones. As 
seen in Table ITT1 WC improves ao significantly over LSD 
and PBE, even much better than TPSS. Note that lattice 
constants should be extrapolated to K to compare with 
DFT results due to thermal expansion. An interesting 
example is ao of cubic PbTi0 3 , which is 3.969 A at 766 K 
(ferroelectric phase transition temperature). It reduces 



TABLE II: Errors (mean absolute relative error) of calculatec 
equilibrium lattice constants ao, bulk moduli Bo, and cohesiv 
energies E c of 18 tested solids, Li, Na, K, Al, C, Si, SiC, Ge 
GaAs, NaCl, NaF, LiCl, LiF, MgO, Cu, Rh, Pd, Ag, at K 
comparing with experiments flOf) - TPSS and PKZB data ar> 
from Ref. pjj. 





LSD 

a 


PBE 


wc 

0-0 


TPSS 

a 


PKZB 

a 


Error (%) 


1.74 


1.30 


0.29 


0.83 


1.65 
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7.6 


8.0 
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Error (%) 


15.2 


5.1 


5.2 







TABLE III: Equilibrium structural parameters for two ferro 
electrics: tetragonal P4mm PbTiC>3 (PT) and rhombohedra 
R3m BaTiOs (BT). The atom positions u z are given in term 
of the lattice constants. 







LSD 


PBE 


WC 


Expt. 


PT 


Vo (A 3 ) 


60.37 


70.54 


63.47 


63.09" 




c/a 


1.046 


1.239 


1.078 


1.071" 




u z (Pb) 


0.0000 


0.0000 


0.0000 


o.ooo' 




u,(Ti) 


0.5235 


0.5532 


0.5324 


0.538 6 




w^OiOi) 


0.5886 


0.6615 


0.6106 


0.612'' 




«*(0 3 ) 


0.0823 


0.1884 


0.1083 


0.112 6 


BT 


v (A a ) 


61.59 


67.47 


64.04 


64.04 c 




a 


89.91° 


89.65° c 


89.86° 


89.87°' 




u z (Ba) 


0.0000 


0.0000 


0.0000 


0.000° 




u,(Ti) 


0.4901 


0.4845 


0.4883 


0.487° 




Wz(OiOi) 


0.5092 


0.5172 


0.5116 


0.511° 




w z (Q 3 ) 


0.0150 


0.0295 


0.0184 


0.018° 



"Low temperature data, Ref. l24l 
'Room temperature data, Re f. |2i 
Tow temperature data, Ref. 1271 



to 3.93 A at K by extrapolation [Ij. a^ BE = 3.971 
A is 1% larger, whereas aff c = 3.933 A, in excellent 
agreement with the extrapolated data. Also note that 
the zero-point quantum fluctuations are not included in 
DFT calculations, which would expand ao about 0.2%. 
WC also predicts more accurate bulk moduli for these 
materials than LSD and PBE (Table [nj . For cohesive 
energies, WC is nearly as accurate as PBE, and much 
better than LSD ( Table HT|) . These results prove that our 
simple model of GGA is very suitable for solids. 

For the ground-state structures of polarized ferro- 
electrics, the lattice strain must be optimized together 
with atomic positions. Unlike cubic systems, a large vol- 
ume for a polarized ferroelectric material favors large 
strain and atomic displacements, and large strain and 
atomic displacements lead to even larger volumes. This 
causes PBE to overestimate the volume of tetragonal 
PbTiOa by more than 10%, whereas the error is only 3% 
for the cubic structure. Table ITTT1 summarizes the LSD, 
PBE, and WC results of fully relaxed tetragonal PbTiOs 
and rhombohedral BaTiOs. It shows that WC predicts 




FIG. 2: Sublimation energies of Bernal- Fowler periodic model 
of ice Ih. 



highly accurate volumes, strains, and atomic displace- 
ments, whereas LSD and PBE underestimate and over- 
estimate these values, respectively. If their model param- 
eters are fitted to first-principles results using WC, MD 
or MC simulations are expected to determine ferroelec- 
tric phase transition temperatures and other properties 
more accurately 

It is well known that LSD fails to predict the correct 
ground states for certain materials, e.g., magnetic bec 
iron |2^| and a-quartz [2^, and PBE can eliminate this 
error. WC predicts correct ground states for both iron 
and quartz with smaller energy differences than PBE, re- 
sulting in lower transition pressures. For iron, the tran- 
sition (bec to hep) pressure of 10 GPa is rather close to 
experiment, but for quartz (a to stishovite), the WC re- 
sult of 2.6 GPa is not sufficient to correct LSD. Since the 
stishovite phase is much more compact and stiffer than 
the a phase, the phonon contributions to energy could 
increase the energy difference greatly. However by per- 
forming first-principles linear response lattice dynamics 
calculations we find that the vibration zero point energy 
difference is only 0.015 eV/Si02, because the a phase has 
high frequency modes the stishovite phase lacks, in addi- 
tion to low frequency modes. These results indicate that 
a better correlation functional for WC is also needed. 

Another interesting case is the weakly interacting 
molecular bonding systems. We calculated the hydro- 
gen bond strength in ice for a periodic Bernal-Fowler Ih 
model [13,113 The wc sublimation energy of 0.73 eV 
per H2O (the measured data is 0.61 eV [23, excluding 
the zero-point vibration) is much better than LSD of 1.07 
eV, but not as good as PBE of 0.63 eV, as shown in Fig- 
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TABLE IV: Jellium surface exchange energies (ax, in 
erg/cm 2 ) computed using the LSD orbitals and densities. 
LSD and PBE values are from Ref. 0, and TPSS values 
are from Ref. [Toll . The last row is the mean absolute relative 
errors against the exact value |34j|. 



r s (bohr) 


LSD 


PBE 


TPSS 


WC 


exact 


2.00 


3037 


2438 


2553 


2519 


2624 


2.30 


1809 


1395 


1469 


1452 


1521 


2.66 


1051 


770 


817 


809 


854 


3.00 


669 


468 


497 


497 


526 


3.28 


477 


318 


341 


341 


364 


4.00 


222 


128 


141 


141 


157 


5.00 


92 


40 


47 


47 


57 


6.00 


43 


12 


15 


15 


22 


Error (%) 


36.7 


16.7 


9.4 


9.8 
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As argued by Perdew 0, in an extended system such 
as metal surface the exact hole may display a diffuse long- 
tail behavior: an emitted electron's exchange-correlation 
hole can extend significantly back into the interior of the 
metal. As summarized in Table IIVI jellium surface ex- 
change energy ox is severely overestimated by LSD and 
underestimated by PBE, and o^ c are better than both 
cr x SD and ct£ be , identical to cr^ PSS for r s > 3.0 bohr, 
where r s = fx 3 —) 1 / 3 . Because the correlation surface en- 
ergy <7£ PSS is very close to <Tq BB , the mean error of <Jx<? 
against the available most accurate Cxc A+ HH > which is 
1.5%, is comparable to 1.1% of TPSS, better than the 
LSD and PBE errors of ctxc of 2.1% and 4.9%, respec- 
tively. The significant improvement of WC upon PBE 
for the jellium surface energy is because the diffuse ra- 
dial cutoff model better describes metal surfaces than 
the sharp cutoff. The good performance of WC on jel- 
lium surface energies suggests it should also perform well 
on metal vacancies |35j |. 

Finally we compared calculated exchange energies of 
5 noble-gas atoms with the Hartree-Fock (HF) results 
6]. The mean errors for LSD, PBE, RPBE, and WC are 
8.77%, 0.89%, 0.18%, and 2.01%, respectively. Compar- 
ing the magnitude of Fx of these GGAs as illustrated in 
Figure one can conclude that among these choices a 
GGA with bigger Fx predicts better Fx of atoms. Al- 
though WC is constructed for slowly varying densities, 
it improves exchange energies of atoms over LSD signif- 
icantly. Furthermore, a proper treatment of correlation 
functional would make WC XC energies of atoms as ac- 
curate as PBE. 

We have shown that exchange enhancement factors Fx 
constructed from different situations perform with differ- 
ent accuracy for the same system. Using PBE correlation 
functional, WC performs excellently for solids, but it is 
less accurate for atoms than PBE; on the other hand, 
RPBE is excellent for atoms, but poorer for solids than 
PBE. It seems impossible to make a GGA which is more 



accurate than PBE for both solids and atoms simultane- 
ously because Fx is a function only of the reduced gra- 
dient s. Our calculations show that F x also depends on 
the variation of |Vn|. Since high density systems often 
have large variations and low density systems often vary 
slowly, we propose that a GGA having Fx = Fx(s,r s ), 
just like the gradient correlation correction, could be uni- 
versally more accurate for atoms, molecules, and solids 
than PBE. The additional parameters of Fx(s,r s ) can 
be fitted to quantum Monte Carlo simulations for spe- 
cific materials, or determined by constraints at r s — > 
and r s — » oo from other theoretical considerations. In 
addition, a better correlation functional will improve the 
XC energy and potential also, and one can construct a 
correlation more compatible with our functional of ex- 
change than PBE. In this way, the accuracy of a simple 
second rung of the ladder of XC approximations, GGA, 
could approach that of the more complicated third rung 
approximation, meta-GGA. 

We have constructed a new GGA which is more accu- 
rate for solids than any existing GGA and meta-GGA. It 
has a very simple form without any empirical parameters, 
and it is ideal for ab inito calculations of certain mate- 
rials, e.g., ferroelectrics, for which exceptionally high ac- 
curacy is needed. It can be generalized to make a GGA 
more accurate for atoms, molecules, and solids than PBE. 
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